function [model, accval]=variational_EM(data, MAXCOUNT, MAXESTEPITER, MAXMSTEPITER, MaxFun, flag, str)

% @ Ayan Acharya, Date: 4.3.2011
% variational EM for BC3E

model=[];
LOWLIMIT=0.01;
convergence=1000000;
model=init_params(data,flag,str);
countVEM=1;
model.k=data.k;
model.k_hat=data.k_hat;

[value1] = cal_likelihood(model, data);

while (convergence> LOWLIMIT && countVEM<MAXCOUNT)
    
    model    = E_step(model,data, MAXESTEPITER, MaxFun, value1, countVEM);
    [value2] = cal_likelihood(model, data);
    model    = M_step(model,data, MAXMSTEPITER, MaxFun, value2, countVEM);
    [value3] = cal_likelihood(model, data);    
    
    convergence=100*abs((value3-value1)/value1);
    value1=value3;
    value(countVEM)=value1;
    disp('count from V-EM');
    countVEM = countVEM+1
    
end

probassgn  = exp(model.alpha_1)./repmat(sum(exp(model.alpha_1),2),1,model.k);
[garb,ind] = max(probassgn');
accval     = cal_accuracy(model, data);

end
